5.0 Introduction

Map and GIS users are mostly confronted in their work with transformations from one two-dimensional coordinate system to another. This includes the transformation of polar coordinates delivered by the surveyor into Cartesian map coordinates (section 2.5) or the transformation from one 2D Cartesian (x,y) system of a specific map projection into another 2D Cartesian (x,y) system of a defined map projection.

Integration of spatial data into one common coordinate system.

Datum transformations are also important, usually for mapping purposes at large and medium scale. An example, map and GIS users are often collecting spatial data in the field using satellite navigation technology and need to represent this data on published maps on a local horizontal datum.

 

5.1 Changing map projection

Forward and inverse mapping equations (see section 4 on map projections) are generally used to transform data from one map projection to another. The inverse equation of the source projection is used first to transform source projection coordinates (x,y) to geographic coordinates (f,l). Next, the forward equation of the target projection is used to transform the geographic coordinates (f,l) to target projection coordinates (x’,y’). The first equation takes us from a projection A into geographic coordinates. The second takes us from geographic coordinates (f,l) to another map projection B. This is illustrated in the figure below.

The principle of changing from one into another projection using the mapping equations.

 

The coordinate system (map projection) of the input data must be known to use the mapping equations for a projection change. If the coordinate system of the input data is not known it may be possible to use a 2D Cartesian transformation. 2D ground control points (GCPs) or common points are then required to determine the relationship between the unknown and the known coordinate system. The transformation may be conformal, affine, polynomial, or of another type, depending on the geometric differences between the two map projections. This is illustrated in the figure below. Refer to section 5.4 on 2D Cartesian transformations for more details.

The principle of changing from one unknown projection into a known projection using a 2D Cartesian transformation. A number of 2D control points are required to determine the relation between both systems.

 

Two-dimensional Cartesian transformations have a different accuracy compared to the transformations based on projection equations. The latter take into account the Earth curvature. This is especially important in the case of large areas and small scale. However, if the control points are coplanar and the extent of the area is not too large, the 2D Cartesian transformation could yield a better model of coordinate relations than the presumed set of projection equations would do.

 

5.2 Datum transformations

A change of map projection may also include a change of the horizontal datum (also called geodetic datum). This is the case when the source projection is based upon a different horizontal datum than the target projection. If the difference in horizontal datums is ignored, there will be no perfect match between adjacent maps of neighbouring countries or between overlaid maps originating from different projections. It may result in up to several hundred metres difference in the resulting coordinates. Therefore, spatial data with different underlying horizontal datums may need a so-called datum transformation. Datum transformations are transformations from a 3D coordinate system (i.e. horizontal datum) into another 3D coordinate system.

Datum shift between two geodetic datums. Apart from different ellipsoids, the centres or the rotation axes of the ellipsoids do not coincide.

 

Suppose we wish to transform spatial data from the UTM projection to the Dutch RD system, and that the data in the UTM system are related to the European Datum 1950 (ED50), while the Dutch RD system underlies the Amersfoort datum. In this example the change of map projection must be combined with a datum transformation step for a perfect match. This is illustrated in the figure below.

The principle of changing from one into another projection combined with a datum transformation from datum A to datum B.

 

The inverse mapping equation of projection A is used first to take us from the map coordinates (x,y) of projection A to the geographic coordinates (f,l) in datum A. Next, the datum transformation takes us from the geographic coordinates (f,l) in datum A to the geographic coordinates (f,l) in datum B. Finally, the forward equation of projection B is used to take us from the geographic coordinates (f,l) in datum B to the map coordinates (x’,y') of projection B. A height coordinate (h or H) may be added to the geographic coordinates.

The transformation parameters to take us from one datum system to another datum system are estimated on the basis of a set of selected points whose coordinates are known in both datum systems. If the coordinates of these points are not correct - often the case for points measured on a local datum system - the estimated parameters may be inaccurate. As a result the datum transformation will be inaccurate. This is often the case when we transform coordinates from a local horizontal datum to a global geocentric datum. The coordinates in the local horizontal datum may be distorted by several tens of metres because of the inherent inaccuracies of the measurements used in the triangulation network. These inherent inaccuracies are also responsible for another complication: the transformation parameters are not unique. Their estimate will depend on which particular common points are chosen, and they also will depend on the amount of transformation parameters estimated.

Mathematically a datum transformation can be realized by relating the geographic coordinates (f,l,h) of both datum systems directly, or indirectly by relating the geocentric coordinates (x,y,z) of the datums. The first method is discussed in section 5.2.2, the latter in section 5.2.1.

 

5.2.1 Datum transformations via geocentric coordinates

Datum transformations via the geocentric coordinates (x,y,z) are 3D similarity transformations. Essentially, these are transformations between two orthogonal 3D Cartesian spatial reference frames together with some elementary tools from adjustment theory. This is illustrated in the figure below.

The principle of changing from one datum into another datum via the geocentric coordinates.

 

The three most applied methods for a datum transformation via the 3-dimensional geocentric coordinates are:

  1. the geocentric translation,
  2. the Helmert 7-parameter transformations (position vector or coordinate frame), and
  3. the Molodensky-Badekas 10-parameter transformation.

Refer to the OGP Guidance note 7: Coordinate Conversions and Transformations including Formulas (page 105-110) for the formulas and a detailed explanation of the formulas.

i. The geocentric translation relates two datum systems through three translations. The method applies a shift between the centres of the two geocentric coordinate systems. This shift is defined by the parameters DX, DY and DZ (or Xo,Yo and Zo).

An example, the ITRF (X,Y,Z) coordinates of a point in the state of Baden-Württemberg (Germany) are:

X (ITRF) = 4,156,939.96m

Y (ITRF) = 671,428.74m

Z (ITRF) = 4,774,958.21m

The three translation parameters for transforming the point from ITRF to the Potsdam datum (the local datum of Germany with underlying Bessel ellipsoid) are given as follow: DX= -635.00m, DY= -27.00m and DZ= -450.00m. This set of parameters provided by the National Imagery and Mapping Agency (NIMA) has been calculated using common points distributed throughout Germany and based on the ITRF. Applying these geocentric translations to the given point results in the following Potsdam (X,Y,Z) coordinates:

X (Potsdam) = 4,156,939.96 - 635.00 = 4,156,304.96m

Y (Potsdam) = 671,428.74 - 27.00 = 671,401.74m

Z (Potsdam) = 4,774,958.21 - 450.00 = 4,774,508.21m

 

ii. The Helmert 7-parameter transformation relates two datum systems through a rotation, an origin shift and a scale factor. The transformation is expressed with seven parameters: three rotation angles (a,b,g), three origin shifts (DX,DY and DZ) and one scale factor (s). An example, the ITRF (X,Y,Z) coordinates of the given point in the state of Baden-Württemberg are transformed to the Potsdam datum. The 7 parameters for transforming the point from ITRF to the Potsdam datum are given as follow:a=+1.04sec, b=+0.35sec, g=-3.08sec, DX=-581.99m, DY=-105.01m, DZ=-414.00m and s=1-(8.3 * 10 -6)=0.9999917. This set of parameters provided by the federal mapping organization of Germany was calculated using common points distributed throughout Germany. Applying these parameters to the given point results in the following Potsdam (X,Y,Z) coordinates:

X (Potsdam) = 4,156,305.34m

Y (Potsdam) = 671,404.31m

Z (Potsdam) = 4,774,508.25m

Tool for the conversion from ITRF(WGS84) to Potsdam coordinates (external link):

WGS84 to Potsdam converter

Potsdam to WGS84 converter

The two sets of transformed coordinates agree at a level of a few meters. The difference in X is 0.38m, in Y around 2.57m, and in Z around 0.04m. These differences can be explained because of the different set of transformation parameters. In a different country, the agreement could be at the level of centimetres, or tens of meters and this depends primarily on the quality of implementation of the local horizontal datum.

The Helmert 7-parameter transformation is considered to be reversible, i.e. the same parameter values can be used to execute the reverse transformation. Thus, the reverse parameters can be applied, whereby a=-1.04sec, b=-0.35sec, g=+3.08sec, DX=+581.99m, DY=+105.01m, DZ=+414.00m and s=1 + (8.3 * 10 -6)=1.0000083, to transform the given point from Potsdam coordinates to ITRF coordinates.

The Helmert 7-parameter transformation can be either a position vector transformation or a coordinate frame transformation. Both transformations are based on the same definition of translation and scale parameters, but a different definition of the rotation parameters. The coordinate frame transformation assumes that the rotations are applied to the coordinate reference frame, while the position vector transformation (also called Bursa-Wolf transformation) assumes that the rotations are applied to the point's vector (see OGP Guidance note 7 for details). The only difference between the two methods is that the rotation has changed sign. The position vector transformation is commonly used, but is not universally accepted. It is therefore essential to communicate the related coordinate transformation method when exchanging datum transformation parameters.

 

iii. The Molodensky-Badekas 10-parameter transformation relates two datum systems through a rotation, an origin shift and a scale factor. This is the same as for the Helmert transformation methods, but instead of deriving the rotations about the origin of the geocentric coordinate system, they are derived at a location within the points used in the determination of the parameters. Three additional parameters, the coordinates of the rotation point (Xp,Yp,Zp), are then required. The transformation is therefore expressed with 10 parameters: three rotation angles (a,b,g), three origin shifts (DX,DY and DZ), one scale factor (s) and the coordinates of the rotation point (Xp,Yp,Zp) given in the source geocentric coordinate system. Compared to the Helmert transformation, the Molodensky-badekas provides usually a better approximation, but the transformation is not reversible.

 

5.2.2 Datum transformations via geographic coordinates

Datum transformations via the geographic coordinates directly relate the ellipsoidal latitude (f) and longitude (l), and possibly also the ellipsoidal height (h), of both datum systems. This is illustrated in the figure below.

 

The principle of changing from one datum into another datum via the 3D geographic coordinates.

 

The applied methods for a datum transformation via the 3-dimensional geographic coordinates are:

  1. the geographic offsets,
  2. the Molodensky and Abridged Molodensky transformation, and
  3. the multiple regression transformation.

i. The simplest method uses geographic offsets. It relates both datum systems with only two parameters, the difference in the geographic latitude (Df) and the difference in the geographic longitude (Dl). The ellipsoidal height (h) is mostly not included. The method is only used for purposes where low accuracy can be tolerated. The equation is:

f datum B= f datum A + Df

ldatum B = ldatum A + Dl

 

ii. A second method of directly transforming latitude, longitude, and height is the standard Molodensky and Abridged Molodensky transformation (it should not be confused with the Molodensky-Badedas transformation). The standard equations directly relate ellipsoidal latitude and longitude coordinates and ellipsoidal height of two datums by deriving the geographic coordinate offsets. The abridged form is found by dropping any terms that are second order in small parameters, and the ellipsoidal height (h) is ignored. Refer to the OGP Guidance note 7 (page 112) for the abridged versions of the formulas and a detailed explanation. The simplified equation is:

Df = f (Dx, Dy, Dz, Da, Df)

Dl = f (Dx, Dy, Dz, Da, Df)

Dh = f (Dx, Dy, Dz, Da, Df)

where Dx,Dy,Dz are the geocentric translation parameters (centre offset between the two datums), Da is the difference in the semi-major axes of the target and source ellipsoids, Df is the difference in the flattening of the two ellipsoids.

The Molodensky equations are commonly used to relate ellipsoidal latitude and longitude coordinates and ellipsoidal height of a local geodetic datum to those of the WGS84 datum. An example is illustrated in the figure below.

The principle of changing from a local datum (Bahamas) into the WGS84 datum via the 3D geographic coordinates using Molodensky. Given are the centre offset parameters (Dx, Dy, Dz) , the difference in the semi-major axes (Da) and the difference in the flattening (Df).

 

The Molodensky equations are easy to use to transform GPS coordinates in WGS84 to a local map coordinate system (map projection) with it's own ellipsoid and datum. Simply the ellipsoid parameters and the datum shift values (centre offset between the two datums) are required to relate the two datums.

List of datums and datum shift values

The list defines a datum as: datum name = ellipsoid name, shift X, shift Y, shift Z, name geographic area or as: datum name = ellipsoid name (in this case a subsection [datum name] provides the shift). The datum shift values refer to WGS84. Tools for the transformation of coordinates from one 3D geographic coordinate system into another using the Molodensky methods:

  • Molodensky (program that finds f, l, and h for a point in the WGS84 datum system. The input are the major-axis (a), inverse flattening (f) and the f, l, and h of the point in a local datum system). MolodenskyGeneralized  (program that finds f, l, and h for a point in a target datum system. The input are the major-axis (a), inverse flattening (f) and the f, l, and h of the point in a source datum system). The program uses the abridged version of the molodensky transformation.
  • MolodenskyInverse  (program that finds DX, DY and DZ, the centre offset between two datums. The input are the major-axis (a), inverse flattening (f) and 3 points (f, l, and h) in both the WGS84 and the local system).

The problem with the Molodensky transformation is the limited amount of parameters (Dx, Dy, Dz, Da, Df) that it uses. This is quite satisfactory for small areas, but for larger areas such as large countries or continents (e.g. US with NAD27, Europe with ED50) significant errors occur (tens of meters). It requires another method to provide good datum shifts for these original types of large area datums.

 

iii. A third method uses multiple regression equations. A series of best-fit equations provide the local shifts in latitude and longitude as a function of position. It is commonly used to relate ellipsoidal latitude and longitude coordinates of continental size datums to those of the WGS84 datum and involve polynomial expressions in the two ellipsoidal coordinates which go up to degree 9 for the time being. The coefficients (transformation parameters) are determined on the basis of coordinate differences of a set of selected points whose coordinates are known in both datum systems. The main advantage of this method over Molodensky (often implemented in geo-software) is that a better fit over continental size land areas can be achieved. The simplified equation is:

Df = f (f, l, a1, a2, a3, .......)

Dl = f (f, l, b1, b2, b3, ...... )

Dh = f (f, l, c1, c2, c3, .......)

For an example refer to the OGP Guidance note 7 (page 90: Polynomial transformation for Spain).

 

5.3 Conversions from geographic to geocentric coordinates and visa versa

The three geocentric transformations decribed in section 5.2.1 are usually combined with conversions between the geocentric coordinates (x,y,z) and the ellipsoidal latitude (f) and longitude (l) coordinates and height (h) in both datum systems. This is illustrated in the figure below.

The principle of a datum transformation via the geocentric coordinates. The datum transformation is combined with conversions between the 3D geographic coordinates and geocentric coordinates in both datum systems.

 

5.3.1 Geographic to geocentric conversion

The conversion from the latitude and longitude coordinates into the geocentric coordinates is rather straightforward and turns ellipsoidal latitude (f) and longitude (l), and possibly also the ellipsoidal height (h), into X,Y and Z, using 3 direct equations. If the ellipsoidal semi-major axis is a, semi-minor axis b, and inverse flattening 1/f, then

X = (u+ h) cos f cos l

Y = (u+ h) cos f sin l

Z = [(1- e2) u + h] sin f

where u is the prime vertical radius of curvature at latitude f and is equal to u = a /(1 – e2sin2 Φ)0.5,
e is the eccentricity of the ellipsoid where e2 = (a2 – b2) / a2 = 2f – f2.

An example, the ITRF(f,l,h) coordinates of the given point in the state of Baden-Württemberg are:

f (ITRF) = 48° 46' 59.6564" N

l (ITRF) = 9° 10' 30.6113" E

h (ITRF) = 330.397m

Applying these formulas to the given point in the state of Baden-Württemberg (Germany) results in the following ITRF (X,Y,Z) coordinates:

X (ITRF) = 4,156,939.96m

Y (ITRF) = 671,428.74m

Z (ITRF) = 4,774,958.21m

Tool for the conversion from geographic (or geodetic) to geocentric (or 3D Cartesian) coordinates (external link):

Geographic to Geocentric converter

To apply the required conversion we need the ellipsoidal height. Datum transformation programs implemented in GIS software often simplify this transformation: i.e. ellipsoidal heights (h) are taken equal to 0. It is obvious that such a simplification will affect the accuracy of the datum transformation. Instead of setting the ellipsoidal heights (h) to zero, the mean sea level heights (H) may be used as input.

 

5.3.2 Geocentric to geographic conversion

The inverse equations for the reverse conversion are more complicated and require either an iterative calculation of the latitude and ellipsoidal height, or it makes use of approximating equations like those of Bowring. These last have millimetre precision for 'Earth-bound' points, i.e. points that are at most 10 km away from the ellipsoidal surface (any point on the Earth surface). Refer to the OGP Guidance note 7 (page 76) for the formulas and a detailed explanation of the formulas. An example, the Potsdam (X,Y,Z) coordinates of the given point in the state of Baden-Württemberg are:

X (Potsdam) = 4,156,305.34m

Y (Potsdam) = 671,404.31m

Z (Potsdam) = 4,774,508.25m

These Potsdam (X,Y,Z) coordinates were computed from the ITRF (X,Y,Z) coordinates by applying a Helmert datum transformation (section 5.2.1). Applying the inverse equations to the given point results in the following Potsdam (f,l,h) coordinates:

f (Potsdam) =48° 47' 3.2752" N

l (Potsdam) = 9° 10' 34.3870" E

h (Potsdam) = 278.825m

Tool for the conversion from geocentric coordinates (or 3D Cartesian) to geographic (or geodetic) coordinates (external link):

Geocentric to Geographic converter

 

5.4 2D Cartesian coordinate transformations

2D Cartesian coordinate transformations can be used to transform 2D Cartesian coordinates (x,y) from one 2D Cartesian coordinate system to another 2D Cartesian coordinate system. The three primary transformation methods are:

  1. the conformal transformation,
  2. the affine transformation, and
  3. the polynomial transformation.

i. A conformal transformation is a linear (or first-order) transformation and relates two 2D Cartesian coordinate systems through a rotation, a uniform scale change, followed by a translation. The rotation is defined by one rotation angle (a), and the scale change by one scale factor (s). The translation is defined by two origin shift parameters (xo,yo). The equation is:

X' = s X cos(a) - s Y sin(a) + xo

Y' = s X sin(a) + s Y cos(a) + yo

The simplified equation is:

X' = aX - bY + xo

Y' = bX + aY + yo

where a=s cos(a) and b=s sin(a). The transformation parameters (or coefficients) are a,b,xo,yo.

 

ii. An affine transformation is a linear (or first-order) transformation and relates two 2D Cartesian coordinate systems through a rotation, a scale change in x- and y-direction, followed by a translation. The transformation function is expressed with 6 parameters:one rotation angle (a), two scale factors, a scale factor in the x-direction (sx) and a scale factor in the y-direction (sy), and two origin shifts (xo, yo). The equation is:

X' = sx X cos(a) - sy Y sin(a) + xo

Y' = sx X sin(a) + sy Y cos(a) + yo

The simplified equation is:

X' = aX - bY + xo

Y' = cX + dY + yo

where the transformation parameters (or coefficients) are a,b,c,d,xo,yo.

 

The diffference between a conformal and an affine transformation is illustrated in the figure below. Both are linear transformations which means that the lines of the grid remain straight after the transformation.

1

a) The uniform scale change of the conformal transformation retains the shape of the original rectangular grid. b) The different scale in x and y-direction of the affine transformation changes the shape of the original rectangular grid, but the lines of the grid remain straight.

 

iii. A polynomial transformation is a non-linear transformation and relates two 2D Cartesian coordinate systems through a translation, a rotation and a variable scale change. The transformation function can have an infinite number of terms. The equation is:

X' = xo + a1X+ a2Y+ a3XY + a4X2+ a5Y2+ a6X2Y+ a7XY2+ a8X3+........

Y' = yo+ b1X+ b2Y+ b3XY + b4X2+ b5Y2+ b6X2Y+ b7XY2+ b8X3+........

Polynomial transformations are sometimes used to georeference uncorrected satellite imagery or aerial photographs or to match vector data layers that don't fit exactly by stretching or rubber sheeting them over the most accurate data layer. The figure below shows a grid with no uniform scale distortions. It may occur in an aerial photograph, caused by the tilting of the camera and the terrain relief (topography). An approximate correction may be derived through a high-order polynomial transformation. The displacements caused by relief differences can be corrected using a Digital Elevation Model (DEM).

5

A grid with no uniform scale distortions. An approximate correction may be derived through a high order polynomial transformation.

 

2D Cartesian coordinate transformations are generally used to assign map coordinates (x,y) to an uncorrected image or scanned map. The type of transformation (usually an affine transformation) depends on the geometric errors in the data set. This is illustrated in the below figure (a). After georeferencing, the image can be aligned (rectified) so that the pixels are exactly positioned within the map coordinate system (figure (b)). For each image pixel in the new coordinate system, a new pixel value has to be determined by means of an interpolation from surrounding pixels in the old image. This is called image resampling.

6

a) Coordinates are assigned to a raster image by means of a 2D transformation using a set of control points. The image is georeferenced; b) The georeferenced image is rectified to match it with the map coordinate system.

 

The parameters (or coefficients or unknowns) of a conformal, affine or polynomial transformation are usually computed with ground control points (GCPs) or common points (also called tie points) such as corners of houses or road intersections, as long as they have known coordinates in both systems. The number of control points required to determine the 4 parameters (a,b,xo,yo) of a conformal transformation must be at least 2. An affine transformation requires at least 3 control points to determine the 6 parameters (a,b,c,d,xo,yo), and 6 control points are required to determine the 12 parameters (xo,a1-a5,yo,b1-b5) of a simple second-order polynomial transformation..

If there are more control points (at least one) than actually required to determine the parameters, we can compute the Root Mean Squares Error (RMSE) using the least squares adjustment. The RMSE is a statistical measure of accuracy for a coordinate measurement, similar to standard deviation, indicating the spread of the measured values around the ‘true’ value. It is a way to quantify the difference between the computed positions and the true positions (residuals) and it represents the average deviation of all the measured positions from their true positions. It is used as an independent check of the quality of the transformation. If the RMSE is not acceptable, then the transformation can be repeated using a different selection or more control points, or a different transformation method can be used.

 

5.5 Summary

In the previous sections we discussed several types of coordinate transformations: projection change, 3D datum transformations and 2D Cartesian transformations. An overview of these transformations is given in the figure below.

The illustration shows how the inverse and forward mapping (or projection) equations are used to transform coordinates from one map coordinate system (projection A) to another (projection B). The projection change may include a datum transformation from one datum (datum A) to another datum (datum B) in case the source projection (projection A) is based upon a different horizontal datum than the target projection (projection B). The datum transformation may take place via a 3D geocentric transformation or directly via a 3D geographic transformation. Alternatively, 2D Cartesian transformations may be used to transform coordinates from one map coordinate system to another (e.g. in case the projection of the input map coordinates is unknown).

Overview of several types of coordinate transformations.

 

5.6 Main references

Knippers, R.A and Hendrikse J. Coordinate transformations. Kartografisch Tijdschrift, KernKatern 2000-3, 2001.

de By, R.A. (editor), Georgiadou, P.Y., Knippers, R.A., Kraak, M.J., Sun, Y., Weir, M.J.C. and van Westen, C.J. Principles of geographic information systems (Chapter 4.2 on spatial referencing), 2nd edition, ITC Educational Textbook, ITC, Enschede, 2001.